Two interacting particles at the metal-insulator transition 



Andrzej Eilmes 1 , Uwe Grimm 2 , Rudolf A. Romer 2 , and Michael Schreiber 2 
department of Computational Methods in Chemistry, 
Jagiellonian University, Ingardena 3, 30-060 Krakow, Poland 
2 Institut fiir Physik, Technische Universitat, D-09107 Chemnitz, Germany 
(Revision : 1.15; compiled February 1, 2008) 

Abstract 

To investigate the influence of electronic interaction on the metal-insulator 
transition (MIT), we consider the Aubry- Andre (or Harper) model which de- 
scribes a quasiperiodic one-dimensional quantum system of non-interacting 
electrons and exhibits an MIT. For a two-particle system, we study the ef- 
fect of a Hubbard interaction on the transition by means of the transfer- 
matrix method and finite-size scaling. In agreement with previous studies 
we find that the interaction localizes some states in the otherwise metallic 
phase of the system. Nevertheless, the MIT remains unaffected by the inter- 
action. For a long-range interaction, many more states become localized for 
sufficiently large interaction strength and the MIT appears to shift towards 

smaller quasiperiodic potential strength. 
71.30.+h, 71.27.+a, 72.15.Rn, 71.23.Ft 
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I. INTRODUCTION 



The physics of the metal-insulator transition (MIT) continues to be at the center of 
current research activities. For two decades it has been known from the scaling hypothesis 
of localization [|l|] that generically a disorder-driven MIT p[ in a free electron system only 
occurs in more than two spatial dimensions, whereas in one or two dimensions an arbitrarily 
small disorder will localize the electronic wave functions. The relevance of many-particle 
interactions for the MIT is much less understood |||4]] . Here we consider the perhaps simplest 
tractable model of an interacting system at the MIT. Namely, we study the case of just 
two interacting particles (TIP) in a particular one-dimensional (ID) quasiperiodic (QP) 
potential. For a single particle (SP) this QP model exhibits an MIT as a function of the 
non-random QP potential strength. 

The problem of TIP in a ID random potential, where the wave functions are always 
localized such that there is no MIT, has already been studied in much detail |5|-|l3|. It 
was argued that a Hubbard onsite interaction U dramatically reduces the localization of 
TIP pair states in comparison with non-interacting and unpaired particles. In particular, 
Shepelyansky proposed an enhancement of the TIP localization length A2 independent 
of the statistics of the particles and of the sign of the interaction such that 



HU) - £/ 2 | (1) 



in the band center with k = 2. Here, Ai is the SP localization length in ID JTj] and U 
is given in units of the nearest-neighbor hopping strength. Microscopic support for the 
derealization was given afterwards by Frahm et al. ||, who observed a behavior A2 ~ A*' 65 
in a numerical investigation employing the transfer- matrix method (TMM). Other direct 
numerical approaches to the TIP problem have been based on the time evolution of wave 
packets [p|,|i5f|, exact diagonalization |K|, Green function approaches P JT2]fl6[] , and TMM 
fnyiTfl . In these investigations an enhancement of A2 compared to Ai has usually been found, 
but the quantitative results tend to differ both from the analytical prediction ([]]), and from 
each other. 



Two of us JTT|] recently studied the TIP problem by TMM but at larger system sizes M 



than Ref. [§] and found that (i) the enhancement A2/A1 decreases with increasing M, (ii) 
the behavior of A2 for U = is equal to Ai in the limit M —>■ 00 only, and (iii) for U 7^ the 
enhancement A2/A1 also vanishes completely in this limit. This raises serious questions about 
the validity of the TMM approach to TIP and in fact it has been argued very recently [|18|| , 
that the TMM approach of Ref. PlJTll may systematically underestimate the localization 
length of a pair state, since it automatically measures a mixture of localization lengths 
originating also from unpaired states. Thus in this work, we will use the TIP-TMM not as 
a tool to extract information about the pair states only, but rather aim at describing the 
general influence of the presence of one particle onto the transport properties of the other. 
At present, it seems well-established by Green function methods P JT^JT6[| that an en- 



hancement A2 > Ai exists, although the validity of Eq. ([]]) is still under debate: the values of 
the exponent k obtained by numerical methods |5|,|8|- |T^ . fr5| -[T7[| range from 1 to 2. In spite of 
these numerical differences, we nevertheless believe that the TIP approach can give mean- 
ingful insights into the interplay of disorder and interaction [[H]]. In particular, the effects of 
interaction on the disorder-driven Anderson transition should be quite interesting already 
for TIP. However, as mentioned above, the disorder-driven MIT requires more than two 
spatial dimensions and so the numerical efforts are close to being prohibitive when including 
interactions. 

Fortunately, the QP — and thus fully deterministic — Aubry- Andre (AA) model fT9| 



exhibits an MIT even in ID, in dependence on the strength of the quasiperiodic potential. 
This model is closely related to the problem of a SP on a 2D lattice in a magnetic field in 
which context it is also known as the Harper model pO|j . At the MIT, the spectrum exhibits 



the famous Hofstadter butterfly shape |2T|, and the spectral and localization properties 
have been studied in great detail ||22|| . In the mathematical literature, the same model is 
also known and studied as the almost-Mathieu equation |23| . 



For this ID model, we can use the TIP approach in a straightforward way in order 
to investigate the effect of the interaction on the transition. Previous studies based on 



perturbative expansions in U and numerical computations of participation numbers in the 



AA model |23| concluded that interaction can lead to the appearance of localized states 
in the metallic regime for TIP. However, although participation numbers are a useful tool 
for characterizing localization properties of states, they may give ambiguous results: in 
some cases, states which are extended or critical may appear to be more localized and vice 
2q| . Moreover, for interacting particles the generalization of localization criteria like 



versa 



the participation number is not straightforward [f26 ]. Thus in this work we concentrate on 
direct calculations of the TIP localization length in the AA model using the TMM for finite 
system sizes. In addition to the onsite interaction, we also consider a long-range interaction. 
Employing the finite-size-scaling (FSS) approach p7 |, we then construct scaling curves from 
which we deduce the localization properties of the infinite system. We find that within the 
accuracy of our results, the critical behavior is not affected by the interactions. But it seems 
that the long-range interaction shifts the critical QP potential strengths towards smaller 
values, thus giving a tendency towards localization. 

The paper is organized as follows. In section [TT], we define the TIP version of the AA 
model and introduce our notations. Section |TTT| reviews the power-series variant of the TMM, 
and the concepts of FSS. In section [TV] we explain the use of a phase-shift parameter in the 
QP potential for reducing statistical fluctuations in the localization length data. Results 
obtained from FSS of the localization lengths for Hubbard and long-range interactions at 
energy E = are presented in section [V|. In section [VI|, we show the localization properties 



of all states of the spectrum. We summarize and conclude in section VII. 



II. THE TIP VERSION OF THE AUBRY- ANDRE MODEL 

The Schrodinger equation for the SP AA model is given as 

0n+l = (E - Hn)4> n ~ 4>n-l- (2) 

Here <p n is a SP wave function, fi n = 2\i cos(an+ (3) is the QP AA onsite potential of strength 
\i with a/27r an irrational number, which we have chosen as the inverse of the golden mean 
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a/27r = (y/E - l)/2, and (3 is an arbitrary phase shift. We remark that a/2n; may be 
approximated by the ratio of successive Fibonacci numbers 1, 2, 3, 5, 8, 13, .... In Fig. [1] we 
show typical data for the SP localization length Ai obtained by TMM for various system 
sizes given by some of the Fibonacci numbers ||28|| . In agreement with previous studies ||19|| , 
this figure suggests already that the MIT occurs ar \x = 1. Of course, further analysis 
like FSS would be necessary for a comprehensive study of this MIT. Here we note that 
in contradistinction to the MIT in the usual Anderson model with onsite random potential 
disorder, in the AA model all states are either extended (/i < 1), critical (/x = 1), or localized 
(/i > 1), and thus no mobility edge, i.e., no MIT in dependence on energy exists. 

In principle, there are many possibilities to extend the SP Schrodinger equation to TIP. 
In order to be most compatible with the TIP approach of Shepelyansky H , we will consider 
a TIP Hamiltonian with an additional QP onsite potential on a chain of length M given as 

M M M 

H = ^2(c n+1 C n + h.C.) + ^2 ^2 U n ,m c n c n c m c mi 
n=l n=l m=l 

M 

+ VnCic n (3) 

71=1 

where and c n are the creation and annihilation operators for the electron at site n and we 
assume that the TIP have different spins. U nym denotes the interaction between particles: 
U ntm = US nm for Hubbard onsite interaction or U n ^ m = U/(\n — m\ + 1) for long-range 
interaction. 



III. THE TRANSFER-MATRIX APPROACH TO TIP 

The TIP Schrodinger equation reads 

1 - 1pn-l,m, (4) 

with ipn^m a TIP wave function which at U = may be written as a product of SP wave 
functions (f) n and <p m . We can rewrite Eq. (f|) in the TMM form similar to a 2D Anderson 
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model on an M x M lattice as {ipn+i, ipn) T 
matrix 



T n (ip n ,ip n -i) T with the symplectic transfer 



T 



■1 



V 



(5) 
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Et-xn- H ± 
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describing the evolution of the wave vectors for the first (n) particle (corresponding to the 
longitudinal direction in the 2D SP TMM approach). Here ip n = (ip ni i, . . . , ip n ,m, ■ ■ ■ , V^Af) 
is the wave vector of slice n, H± is the SP hopping term for the second (m) particle (cor- 
responding to the transverse direction) and (Xn)i,m — [fJ-n + + U ntm }5 itm codes the QP 
potential and the interaction || . Note that in this approach the symmetry of the wave func- 
tion remains unspecified and we cannot distinguish between boson and fermion statistics. 

The evolution of the state is determined by the matrix product = E[n=i T n and we 
have 

/ lh*r., \ ( llu \ 

(6) 

Usually, one studies a quasi-lD system of size M x N with M N . However, in the 
present problem, both directions are restricted to n,m < M and iterating Eq. (^j) only 
N = M times will not give convergence. Frahm et al. j§| have solved this problem in their 
TMM study by exploiting the Hermiticity of the product matrix Qm = Tm^m'- Continuing 
the iteration (||) with r M , then with tm, and so on, until convergence is achieved, yields the 
eigenvalues exp[— 2M^j\ of Qm- This is the well-known power method for the diagonalization 
of Hermitian matrices [p9|| . The smallest positive Lyapunov exponent 7 m i n determines the 
slowest possible decay of the wave function and thus the largest localization length A max = 
1 / 7 min for given energy E and phase shift (3. We now define the localization length A of the 
two-particle wave function ij) n>m as A max of the transfer matrix problem (|) and expect it to 
reflect the influence of the particle interaction. 

According to the one-parameter scaling hypothesis [I]] , which has been verified with very 
high accuracy for random potentials \i n |2"7|] , the reduced localization lengths A(M)/M scale 
onto a single scaling curve, i.e., 
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\(M)/M = f(£/M). (7) 

For the AA model considered here, we are not aware of any previous FSS study. Indeed, it 
is not a priori obvious that one-parameter FSS should be valid for the AA model. At least 
for a given single phase shift f3, it is clear from Fig. [1| that we need to go to rather large 
system sizes in order to suppress the fluctuations around /i = 1 and to be able to use the 
FSS approach. However, as we will explain in the next section, we may use different values 
of (3 as being analogous to the different disorder realizations in the Anderson model. As 
usual, we may then determine the finite-size-scaling (FSS) function / and the values of the 



scaling parameter £ by a least-squares fit |27 . 



IV. AVERAGING OVER DIFFERENT (3 

The localization length calculated for given system size and QP potential // depends 
significantly on the (3 value as shown for SP in Fig. |2|. This means that the decay length 
varies depending on the phase shift of the potential along the chain. One may expect that 
the chain length M will also influence the results by changing relative phases of the potential 
at the ends. 

Therefore we have restricted our calculations to the chain lengths given by the Fibonacci 
numbers mentioned in section [□], because for our choice of a, this assures that the phase 
difference of the potential at both ends of the chain will be similar, i.e., approaching zero 
with increasing M. We note that our numerical results presented in the next section do not 
change significantly, when we alternatively use the rational approximants for a/27r instead 
of the irrational number defined in section 0. 

Still, the dependence of Ai on the system size M for a given value of (3 shows much 
structure which makes simple extrapolations towards the infinite system or FSS impossible. 
This dependence is also responsible for fluctuations of the SP Ai close to fi = 1 which are 
visible as peaks in Fig. [I] for small Fibonacci number M [2S|. Only for very large M, the 



fluctuations become small. On the other hand, finite systems with different values of (3 may 



be viewed as different parts cut out of the infinite QP model. This then suggests that we may 
reduce the fluctuation effects by averaging over many such small pieces or, equivalently, many 
different values of randomly chosen (3. Thus different (3 values are analogous to different 
disorder configurations used in the Anderson Hamiltonian. Fig. |3] presents such an average 
over 1000 /3 values for the SP localization length. As expected, the fluctuations visible in 
Fig. [I] disappear even for small systems and extrapolations to large M and FSS are now 
possible. 



V. LOCALIZATION PROPERTIES AT E = 

We now turn our attention to the problem of TIP and study the effects of interaction 
on the localization lengths obtained by TMM, restricting ourselves to E = 0. To this end, 
we have computed the localization lengths for 6 system sizes M — 8, 13, 21, 34, 55, and 89, 
for 80 QP potential strengths /x ranging from 0.56 to 4, and for 6 interaction strengths 
U = 0, 0.5, 1, 1.5, 2, and 10, with onsite and with long-range interaction. Typically, for each 
such triplet of parameters (M, /i, U) we averaged over at least 1000 different (3 realizations. 
We note that as for the case of TIP in a random potential [ftp ], attractive and repulsive 



interaction strengths give the same results at E = and we can thus restrict ourselves to 
U > here. 



A. Hubbard interaction 



Figs. H and |5] show the FSS results for /3-averaged data at energy E = for onsite 
interaction strength U = and U = 1. As can be seen, the coalescence of data for various 
values of /x is not perfect and in fact certainly worse than, e.g., for onsite random disorder |2"? 



This is especially visible on the extended side /x < 1. Nevertheless, the figures clearly show 
the existence of two branches of the scaling curve as in the 3D Anderson model []2?|]. This 
indicates, in agreement with the above considerations for the SP AA model, the presence 
of localized states for /x > 1 and extended states for /x < 1. The MIT appears at a critical 
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QP potential fi c which is close to 1. The values determined from the FSS procedure are 
/i c = 1.01 ± 0.02 for U = and /i c = 1.04 ± 0.04 for £7 = 1. Corresponding FSS plots for 
U = 0.5, 1.5,2 and 10 all consistently give \x c ~ 1. We attribute the small deviations from 
the critical value fi c = 1 of the SP case to the fluctuations in the data. 

Thus the MIT does not get shifted by the Hubbard interaction and the transport prop- 
erties of one particle in the presence of another remain unchanged. On the metallic side 
of the transition (// < 1) this is immediately clear: the interaction is supposed to localize 
0(M) TIP states out of the 0(M 2 ) states in the unsymmetrized Hilbert space |24]. The 



TMM inherently measures the longest localization length and thus simply misses the few 
shorter localization lengths induced by the interaction. On the localized side, however, we 
could expect the interaction to delocalize these TIP states which might be visible even by 
TMM. However, as discussed in section |, this effect is not present in the TIP-TMM or at 
least too small to be visible ||18|| . 

The scaling parameters £ obtained by FSS according to Eq. ([7|) are expected to diverge 
at the transition as £ ~ |/x — n c \~ v with the critical exponent v. In Fig. [(] we show the 
dependence of £ on the QP strength \x. The divergence at /i c ~ 1 is clearly visible. A power- 
law fit gives v = 0.8 ± 0.2 both for U = and U — 1. The large error of the estimate is 



due to the fluctuations in the data near the critical point p7| . Furthermore, in the localized 



regime of the SP AA model is has been shown that JT9| 



Ai~l/]n(l + |/x-jx c |), (8) 

which yields v = 1 by expansion around /i c . In order to check whether this equation 
holds also for TIP we examined the dependence of l/£ on In (1 + |/i — fi c \). The results are 
displayed in Fig. 0. The slope of the best fit line is 1.00 ± 0.02 for U = and 1.01 ± 0.03 
for U — 1. This suggests that onsite interaction does not change the critical behavior at the 
MIT. 
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B. Long-range interaction 



We now consider the long-range interaction denned in section O. The FSS plot for U = 1 
in Fig. |8| is qualitatively the same as for Hubbard interaction. We find localized states for 
fi ^> 1 and extended states for fi <C 1. In Fig. |], we have included the variation of the scaling 
parameter £ with \x for this case. The divergence of £ occurs at fi c = 0.92 ± 0.04 indicating 
that the MIT has been shifted towards smaller values of the QP potential strength fi. FSS 
plots for U = 0.5, 1.5 and 2 suggest that this shift becomes somewhat more pronounced for 
larger U, decreasing to \i c ~ 0.9 for U = 2. 

This behavior may be rationalized by keeping in mind that for a long-range interaction, 
contrary to the case of Hubbard interaction, all states will eventually feel the interaction- 
induced tendency towards localization on the extended side of the MIT, as we will show 
for small systems in the next section. Thus even the most delocalized states at E — will 
become more localized for sufficiently large U. However, in order to answer the question 
whether long-range interaction indeed shifts the MIT towards weaker QP potential strength, 
additional calculations with still higher accuracy would be necessary. These require, however, 
a prohibitive numerical effort when using the present power-series method. 

The critical exponent for U — 1 calculated as in section |V A| is v = 1.0 ± 0.2 and the 
respective slope in Fig. [7| is 0.97 ± 0.03. These values are compatible with our results for 
onsite interaction within the error limits. Therefore, the critical behavior is similar to the 
SP case and onsite interaction. 

VI. ENERGY DEPENDENCE OF THE LOCALIZATION PROPERTIES 

In the previous section, we have rationalized the persistence of the MIT in the presence 
of interactions by assuming that on the extended side the onsite interaction localizes a small 
number of states leaving the rest unaffected. To further examine this effect with TMM we 
calculate the dependence of A on the energy E for a single value of (3 and a small system 
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size. 



A. Hubbard interaction 

Fig. ^| presents results for the inverse localization length A -1 obtained by TMM on the 
metallic side. Also shown are the values of the eigenenergies Ei. The TMM accurately 
shows that transport at energies not corresponding to eigenstates is suppressed, because the 
incoming wave function decays exponentially. On the other hand, A _1 (-E) decreases rapidly 
towards zero when E is approaching an eigenvalue Ei as shown in Fig. ||. This has also been 
observed in the SP case For U — 0, we find a few cases where A -1 remains large even 
at the energy of an eigenstate. From an analysis of the corresponding wave functions, we 
can identify these states with boundary states where the particles are localized close to the 
ends of the finite chains. 

The comparison of the plots for U — and U = 1 shows that, while the energy of most 
states changes only slightly, there are a few states which move to significantly larger energies. 
Their localization lengths are apparently much shorter. The calculation of /^-averaged decay 
lengths for different M shows that states at the verge of the spectrum at E = 4.6 remain 
extended for U = 1 while the states at E = 5.3 are localized. There are also some states 
within the main part of the spectrum which shift to higher energies. Some of them are 
visible in Fig. ^ as they enter the energy gaps. For sufficiently strong interaction U — 8 
there are 13 localized states which split off the remaining spectrum. The calculations for 
other system sizes support the conclusion that the interaction localizes M out of M 2 states 
for system size M. These states correspond to both particles residing on the same site and 
interacting via the Hubbard U. The other states remain extended and do not change their 
energy significantly. 

In the localized regime (/i > 1) the interaction has a similar effect, i.e., for sufficiently 
large U it shifts M states above the main part of the spectrum and increases their localiza- 
tion; the remaining unshifted states also stay localized. These results are in agreement with 
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those of Ref. fM \, when we keep in mind that our numerical method does not allow us to 
see accurately an eventual derealization at intermediate U. 



B. Long-range interaction 

Fig. |It] presents respective results obtained for long-range interaction. Again, the in- 
teraction shifts states to higher energies and shortens their decay lengths. However, as the 
particles feel the interaction at any separation, all states change their energy in agreement 



with section |VB| . This is especially pronounced in Fig. [10| for U = 10. The most prominent 
shift is the change at the high energy part of the spectrum. For extremely large interaction, 
e.g. U = 1000, the spectrum splits into M groups of states reflecting the number of sites 
at which two particles may reside at given separation, i.e., for system size M there are M 
states for separation n — m = 0, and 2M — 2 states for separation \n — m\ = 1, and so on. 

VII. CONCLUSIONS 

We have demonstrated that it is possible to perform FSS for a system with two interacting 
particles in a ID QP Aubry-Andre (or Harper) potential. We find two branches in the FSS 
curves which correspond to localized and extended behavior. The roughness of the FSS plot 
is probably an effect of small system sizes and insufficient averaging and should disappear 
for larger systems, requiring, however, much larger computational effort. On the other hand 
it may be that one-parameter scaling is not strictly valid in this QP model as evidenced by 



the results for the even chain lengths M = 34 and 144 p8[ . Nevertheless, even in this case 
the presence of localized and extended branches as in Figs. f|, || § indicates the existence of 
an MIT. 

The FSS results for energy E = show that the MIT exists in these TIP systems both for 
the non-interacting and the interacting case. The transition point /i c does not depend on the 
Hubbard interaction strength U and is located at QP potential strength fi c « 1. However, 
a large enough long-range interaction shifts fi c towards smaller QP strength. Within the 
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numerical accuracy of our data, the critical behaviour of the localization length is not affected 
by the Hubbard and the long-range interaction. 

The dependence of the decay length on the energy as calculated by TMM confirms the 
results obtained by other methods [[OJ that a large enough interaction localizes pair states 
simultaneously increasing their energy and leaves the rest of the states almost unaffected. 

In closing we remark that our results may also be viewed independently of the TIP 
problem, by noting that the present problem of two particles in a ID QP potential may also 
be seen as SP problem in a particular realization of a 2D QP potential. Similar systems 



have been investigated previously, e.g., in Ref. 30 within the Landauer approach 
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FIG. 2. Inverse of the localization length Ai for the SP Aubry-Andre model as a function of 
phase shift (3 for E = and M = 13. Different symbols indicate QP potential strength \i = 2 (V), 
1.5 (*), 1 (□), and 0.5 (A). 



18 




0.5 1.0 1.5 2.0 



FIG. 3. Localization length Ai for the SP Aubry-Andre model as a function of QP potential 
strength \i for E = 0, averaged over 1000 /3-values. The system size is increasing from bottom to 
top. Note the MIT at fi = 1. 
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FIG. 6. Scaling parameter £ as a function of QP potential strength fx for U = (o), for onsite 
interaction with U = 1 (□), and for long-range interaction with U = 1 (O). 
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FIG. 7. Inverse scaling parameter l/£ as a function of QP potential strength \i as in Eq. (|8|) 
for U = (o), for Hubbard interaction with U = 1 (□), and for long-range interaction with U = 1 
(O), consecutively shifted by 1 for clarity. The lines indicate linear regression fits to the data in 
the localized regime. 
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FIG. 8. Scaling function (0) for long-range interaction U = 1, E = and various /i. The 
characters are chosen as in Fig. 
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FIG. 9. Inverse localization length as a function of energy for a QP potential strength fi = 0.9 
at (3 = y/2, M = 13 for U = and for two Hubbard interaction strengths U. Plots for different U 
are vertically shifted for clarity. The eigenenergies are indicated by (+). 
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FIG. 10. Inverse localization length as a function of energy for a QP potential strength [i = 0.9 
at (3 = y/2, M = 13 for U = and two long-range interaction strengths U. Plots for different U 
are vertically shifted for clarity. The eigenenergies are indicated by (+). 
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